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I. INTRODUCTION 

For more than two decades the correlation induced 
metal-insulator transition (MIT) and its characteristics 
has been one of the challenging problems in condensed 
matter physics This MIT is often accompanied by a 
symmetry breaking and the development of long range 
order In dimension D=l this ordering can only 

be connected with the breaking of a discrete symmetry. 
Examples include commensurate charge density waves 
(CDWs) and Peierls dimerization phenomena. The MIT 
can be driven either by varying the electron density or 
the electron-electron interaction strength. If the MIT is 
approached from the metallic side, the electrical conduc- 
tivity or the electronic compressibility may be used to 
characterize the transition; the approach of the transi- 
tion from the insulating side may on the other hand be 
properly captured by the divergent behavior of the elec- 
tric susceptibility |3|. 

In addition, the features of the insulating phase may 
change qualitatively with the interaction strength, al- 
lowing for quantum phase transitions between insulating 
phases. The extended Hubbard model at half-filling with 
an on-site (C/) and a nearest neighbor (V) Coulomb re- 
pulsion provides a prominent example with a transition 
from a Mott insulator (MI) to a CDW insulator in the 
vicinity of the U — 2V line in the phase diagram Q . Re- 
markably, the transition involves an intermediate phase 
with a bond-order wave (BOW) |,|. 

In recent years particular attention has been given to 
another example for an extension of the Hubbard model 
which includes a staggered potential term . The cor- 
responding Hamiltonian has been named the "ionic Hub- 
bard model" (IHM); at half-filling this model undergoes 
a different kind of insulator-insulator transition from a 
band to a correlated Mott-Hubbard type insulator. Con- 
flicting results have so far been reported regarding the 
nature of the transition, the possibility of two rather 
than one critical point, or the appearance of BOW or- 
der [To|J^,P,pd]-p^ . Given the numerous unresolved issues 



we reinvestigate in detail the IHM using both, numeri- 
cal and analytical tools. Specifically we study the ID 
Hamiltonian 

H = -t ^(1 + i-iyS) (aLc,+i, + H.c) 

i ia 

where creates an electron on site i with spin cr, 
rij^ = cjg.Cjj^. A is the potential energy difference be- 
tween neighboring sites, and S a Peierls modulation of 
the hopping amplitude t. In the limit A = 5 — 0, Eq. (|l|) 
reduces to the ordinary Hubbard model, the limit A = 
and 5 > is called the Peierls- Hubbard model (PHM), 
and the limit A > and 5 = is usually referred to as 
the IHM. In the following, we will focus mainly on the ef- 
fect of the on-site modulation A, so we implicitly assume 
d — except where stated otherwise. 

The IHM was first proposed and discussed almost 20 
years ago in the context of organic mixed-stack charge- 
transfer crystals with alternatin g do nor (D) and acceptor 
(A) molecules (. . . D+pA-^p . . .) These stacks form 

quasi-lD insulating chains, and are classified into two 
categories depending on the amount of charge transfer 
p: quasi-neutral for p < 0.5, and quasi-ionic for p > 0.5. 
Torrance et al. p4| found that at room temperature and 
ambient pressure, these materials are either mostly ionic 
(p w 1) or mostly neutral (p w 0), but several systems un- 
dergo a reversible neutral to ionic phase transition (NIT) 
jisf , i.e. a discontinuous jump in the ionicity p upon 
changing temperature or pressure with the simultaneous 
appearance of a Peierls lattice distortion The 
original modelling of the NIT in charge-transfer salts in- 
cluded different Coulomb repulsion parameters on A and 
D molecules as well as intersite interactions [^8[^ and 
electron-molecular vibration coupling [ p0[ . In these ex- 
tended ionic Hamiltonians the NIT was found to be first 
or second order depending on the relative strengths of 
the different interactions. 
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In addition, the IHM has been used in a quite dif- 
ferent context to describe the ferroelectric transition in 
nerovskite materials such as BaTiOs [|l]j2|Jl^ or KNbOa 
p4| . The IHM here is used as a simple model for a bi- 
nary compound, AB, with the A sites representing the 
transition metal element (e.g. Ti), and the B sites repre- 
senting the oxygen. The staggered potential in this case 
is the energy difference of the levels on the A and B sites 

The results of our work on the IHM support the exis- 
tence of a single transition at a critical coupling [4 (A) 
from a band- to a correlated insulator (CI) phase. The 
transition originates from a ground-state level crossing 
with a change of the site-parity eigenvalue. From our re- 
sults we identify the strong coupling phase at [/ > [4 (A) 
as a CI with long-range CDW order. Two possibilities 
arise for the unconventional CI phase: it either has a 
vanishing spin excitation gap and exhibits on top of a 
CDW pattern an identical power law decay of SDW and 
dimer-dimer correlations; or the spin gap remains finite 
in the CI phase with true long range BOW order coex- 
isting with a CDW. Within the achievable numerical ac- 
curacy our data support the latter scenario. Furthermore, 
we provide the first numerical evaluation of the optical 
conductivity for the IHM which allows for an additional 
characterization of the "metallic" transition point [/c(A) 
where the optical gap closes. The distinction between 
charge and optical gaps is crucial for the structure of the 
ground-state phase diagram. 

The outline of our paper is as follows: In sections II and 
III we summarize previously proposed scenarios for the 
ground-state phase diagram using bosonization and sym- 
metry arguments. Section IV contains a detailed reanal- 
ysis of the bosonization approach to the IHM. In sections 
V and VI we present our numerical Lanczos and DMRG 
results, from which we draw conclusions for the ground 
state phase diagram, and in section VII wc discuss dy- 
namical DMRG results for the optical conductivity. Fi- 
nally, we conclude and summarize in section VIII. 



II. INSULATOR-INSULATOR TRANSITION 

A good starting point for understanding the existence 
of a phase transition in the IHM is the atomic limit |^,^ . 
For i = 0, it is immediately seen that at half- filling and 
[/ < A the ground-state of the IHM has two electrons 
on the odd sites, and no electrons on the even sites cor- 
responding to CDW ordering with maximum amplitude. 
On the other hand, for C/ > A each site is occupied by 
one electron. Obviously, for t — a transition occurs at 
a critical value Uc — A. This transition is expected to 
persist for finite hopping amplitudes t > 0, where the al- 
ternating potential still defines two sublattices, doubling 
the unit cell and opening up a band gap A for U = 
at fc = ±7r/2. For t > the critical coupling shifts to 
Uc{t) > A, where Uc increases monotonically with in- 



creasing A. For U,A^t the system is close to the 
atomic limit, and Uc approaches A from above. 

For [7 = the ground state at half-filling is a CDW 
band insulator (BI) , whose elementary spectrum consists 
of particle-hole excitations over the band gap. We con- 
sider a system to be in a BI phase when the criterion 
As = Ac holds, where the spin (A5) and the charge 
gap (Ac) are given by 



As^EoiN 
Ac = EoiN 
+ EoiN 



L 
L 



1) 



-Eo{N = L,S, 
1/2) 

l/2)-2Eo{N 



0), 



L,S, = 0). 



Eo{N,Sz) is the ground-state energy, L the system 
length, N the number of electrons, and Sz the z- 
component of the total spin. For U > Uc and U » i, A, 
the charge gap is set by the Coulomb interaction U, and 
the system is a CI with Ac > Ag. However, in contrast 
to the cases with A = or < = 0, CDW order is expected 
for all finite values of U. 

The IHM is distinctly different from the PHM. This 
is also a BI at [/ = 0, but in contrast to the IHM has 
Ac > As > for any value of [/ > 0, i.e. the phase 
transition from the Peierls BI to the CI occurs at ?7c = 0. 

A renewal of interest in the IHM started with the 
bosonization analysis of Fabrizio et al. (FGN) |^ and 
exact diagonalization studies where a new scenario 
for the intermediate region L/ ~ A was proposed. 
Several subsequent numerical treatments of the model 
p7| , p8 11-1^ led to partially contradicting results. Dif- 
ferent conclusions were reached about the nature of the 
insulator-insulator phase transition, the possible exis- 
tence of two transitions, and the question whether the 
spin gap closes in the strong coupling phase. 

FGN argued that for small but finite U the BI persists 
up to a critical value Uci , where Ac = and the system 
might be "metallic" . Upon further increasing U they pre- 
dicted a "spontaneously dimerized" , (SDI or equivalently 
a BOW) intermediate phase, which should undergo a con- 
tinuous transition into the MI phase at a second critical 
value Uc2 > Uci- A BOW ground state is characterized 
by long-range staggered bond-density correlations 



9B{r) = ^ X! {i^o\BtBi+r \ipo) , 

i 



H.C. 



(2) 
(3) 



and implies a finite spin gap As > 0. 

Evidence for such a BOW state has been reported by 
Wilkens and Martin in a quantum Monte Carlo (QMC) 
study ijllj. More precisely, the authors of Ref. jll] ob- 
served a single transition from the BI to the correlated 
BOW phase. On the contrary, DMRG results ||2|l and 
Lanczos studies of level crossings in the excitation spec- 
tra support the existence of two transitions with an inter- 
vening BOW state p3[. From finite size extrapolations 
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it was furthermore concluded that A5 = above the 
critical region jl^j2^,^. Whether CDW order persists 
also in the strong coupling regime remained unsettled in 
previous work. 



0. For the IHM with > i, A the effective Heisenberg 
spin model 



III. SYMMETRY ANALYSIS 

Insight into the nature of the BI-MI transition is ob- 
tained from a symmetry analysis. The IHM is invariant 
with respect to inversion at a site and translation by two 
lattice sites. Thus, any nondegenerate eigenstate of H is 
also an eigenstate of the operators that generate the cor- 
responding symmetry transformation. If we denote the 
site inversion operator by P, defined through 



H.P^=A-^a fori = 0,...,L-l 



(4) 



and Tj for a translation by j sites, then any nondegener- 
ate eigenstate IV'™) of H must obey P\tpn) = ±1^"") and 
'Z^IV'ti) = IV'n)- Because [H,Ti] 7^ 0, a non-degenerate 
eigenstate of H can not be an eigenstate of Ti. 

For the half-filled Hubbard model (A = 5 = 0) the site- 
parity eigenvalue of the ground state can be determined 
in the limits [/ = and f7 > i (|]. For [/ = 0, the ground 
state is a direct product of spin up and spin down Slater 
determinants, both formed from the same occupied spa- 
tial wavefunctions with the same parity P„ = ±1, so the 
parity eigenvalue of the total wavefunction is given by 
their product P = x = ±1. Hence, the state at 
?7 = is even under site inversion. On the other hand, in 
the large U limit the mapping to the Heisenberg Hamil- 
tonian can be used to show that for L = An with peri- 
odic boundary conditions (PBC) or L = An + 2 with an- 
tiperiodic boundary conditions (APBC) the ground state 
obeys P|V'o) — — IV'o) (for details see In finite chains, 
these combinations of chain lengths and boundary con- 
ditions ensure that /c = 7r/2 is an allowed momentum 
in the Brillouin zone (BZ). This is important at half- 
filling, where k = ±tt/2 is the BZ boundary. It follows 
that upon increasing the number of sites L, the ground 
state for U ^ t will be odd with respect to P as long 
as = 7r/2 is an allowed k value, and this feature will 
persist in the thermodynamic limit L +00. For the 
ordinary Hubbard model the discrete inversion symme- 
try changes therefore between U = and U ^ t. Since 
the model has no phase transition for U > 0, the ground 
state has P = +1 only for U = 0, and P = — 1 for any 
U >0. 

However, in the IHM the phase transition from a BI 
to a CI occurs at some finite Uc > 0. This suggests 
that the parity of the ground state remains even not only 
for U = 0, but for all U < Uc- At Uc, a, ground-state 
level crossing occurs on finite chains, as confirmed by 
exact diagonalization studies (see below) , connected with 
a site-parity change. 

In the strong-coupling limit, the IHM was argued in 
some previous treatments P,H,f3i to be a MI with Ag — 



was derived to describe the low-energy physics 
(E) the exchange couplings are given by 



In Eq. 



J 



4t2 r 
77 



1 



4t^ 1 + 4a;^ 



it^ {1 + 



C/2 (l-x2)3 



U^ (l-x2)3 



(6) 



where x = A/U. It is important to recall that in the 
derivation of the effective spin model the Hilbert space of 
the IHM is divided into two subspaces - one containing 
no doubly occupied sites, and one containing all states 
with at least one doubly occupied site. The atomic limit 
of the IHM is taken as the unperturbed problem Hq, and 
the kinetic term in (|l|) is treated as a pertubation. The 
ground state of Hq is therefore in the subspace without 
double occupancies. The kinetic term lifts the spin de- 
generacy in this subspace. Introducing the projection 
operator P = Y[i ~ ' effective Hamiltonian 

Hcff = P exp{iS)H eKp(—iS)P is obtained by a canoni- 
cal transformation S and an expansion in powers oi t/U 
around Hq [ ^J30| . 

The result (|5|) would imply that in the strong-coupling 
limit of the IHM the low-energy physics is qualitatively 
similar to that of the Hubbard model, with modified ex- 
change coupling constants J and J' . For next-nearest 
neighbor couplings J' < 0.24 J the spin gap vanishes 
thus the system would be a true MI. The coupling con- 
stants (^) satisfy the condition J' < 0.24J at least for 
U > 3M for A < t and U > 3.6A for A > <. However, 
the effective model (^) is invariant with respect to trans- 
lations by one lattice spacing, whereas the original IHM 
is invariant only with respect to translations by two lat- 
tice spacings. Thus, the effective spin Hamiltonian has 
a higher symmetry than the original model from which 
it was derived, and the arguments of the strong-coupling 
expansion in favor of Ag = for U > Uc cannot be 
considered rigorous as was also pointed out in Ref. jll] . 

Of course, also in the Hubbard chain for any U < 00 
the double occupancy is not zero. However, these "vir- 
tual doublons" are distributed equally over all sites, and 
therefore do not change the translational symmetry. If 
only the spin excitations are considered, these doublons 
are never traced. 

As originally discussed in Ref. the existence of a 
level crossing on finite chains at a critical Uc leaves the 
following possibilities for the ground state of the IHM 
at U > Uc in the thermodynamic limit: 1. The ground 
state remains unique with site-parity quantum number 
P = — I. This necessarily implies the absence of BOW 
and the system is in the MI phase with a vanishing spin 
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gap and absence of long range CDW. Remarkably and 
left uncommented before, this scenario implies the highly 
unusual situation that the strong coupling ground state 
has a higher symmetry than the IHM Hamiltonian. 2. 
The ground state is twofold degenerate and supports a 
long-range BOW+CDW order. This possibility requires 
the spin gap to remain finite These intimately con- 
nected properties are violated in the conclusions reached 
in Ref. [^. This second scenario though still allows for 
an additional second transition into the MI phase as char- 
acterized in 1. 



IV. BOSONIZATION RESULTS 

The standard bosonization procedure allows to ex- 
press the Hamiltonian density in the following way 
H = He + Ti-s + Hcs , where the charge and spin degrees 
of freedom are represented by a sine- Gordon model for 
scalar fields (jjc and 0s, respectively: 



+ j dx{|[p,2(x) + (a.</),n} 

J da; sin (^-^/2TrK^(l)c^ cos (^^/2nKs 



A 

TTOo 



i.e. the half-filled Tomonaga-Luttinger (TL) model with 
ionic distortion. The non-interacting system corresponds 
to the particular limit Kc = = 1. For generality we 
do not assume SU(2) symmetry of the spin channel and 
therefore the parameter Kg is not fixed to unity. 

If we decouple the interaction term in a mean-field 
manner by introducing 



= A 



(11) 
(12) 



u 



27raQ 



cos 



U 
27raQ 



(7) 



(8) 



rris = A - (sin(v/27ri^c0c)) , 
the bosonized Hamiltonian reads Ti. = Tic + Tig where 

-^sin(V2^<^c)|, (13) 

-— cos(V2vri^s0s)|- (14) 



ao is the lattice constant. The ionic term determines the 
spin-charge coupling 



Trao 

We can estimate the renormalized masses as 



(15) 
(16) 



He 



y^27r Kc(l>c{x)) cos (V2Trcl)g{x)\ . (9) where A is a c ut-off e nergy. Simila rly, th e expectation 



Here Pc(s)(x) is the momentum conjugate to 4>c(s){x), 
Vc{s) is the velocity of charge (spin) excitations and Kc — 
1 — U /Airt. This means that the model (|l|) is in a strong- 
coupling regime for arbitrary U > 0, and at A = the 
dynamically generated mass determines the charge gap 
in the system. 



A. Forward scattering model with ionic and/or 
Peierls distortion 



We first consider the limit of a weak Hubbard interac- 
tion U <C A <C where the properties of the system are 
mostly determined by the ionic distortion. In the infrared 
(low energy - large distance) limit the Umklapp and back- 
ward scattering processes, described in Eqs. (0)-(|l) by 
the terms t/cos [y/SirKcct'cix)) and t/cos (V87r</'s(a;)) , 
respectively, are frozen out. Therefore, we initially ne- 
glect these terms and consider the Hamiltonian 



values of sin(v27rK7</'c) and cos{y/2'KKs(t>s)) are esti- 
mated as 



(sin(v/2^0,)) c (^) 

(cos(V2^'/>c)) ^ (^) 
The self-consistency equation 



A = 



(17) 

(18) 
(19) 



(sin(V2ii?7</'s)> (sin(\/2^0,)) 
leads to 

X - (x) • (^0^ 

Using Eqs. ([rT|)-([l9|) one easily finds that 



M, = = A 



2/(4-Kc--R's) 



Htl = 



dx{^-^[P^{x) + {d.cf,cf]} 



(10) 



and the system is a BI with a weakly renormalized 
value of the excitation gap. In particular, for the non- 
interacting system Kc — Kg = 1 and Ms = Ale = A. 
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For the PHM the bosonized expression for the Peierls 
distortion term is given by 



J % J % 

+ f / ^[1 - l^cf - cos(2</,(x)) + ...] , (28) 

where dots denote strongly irrelevant terms of higher crit- 
ical dimensionality, wc obtain the following renormaliza- 
tion of the model parameters: 



_S_ 

TTOo 



cos (\/2^(?!)c) cos (y2^(l)s^ . (21) 



By following again the same steps as above the half-filled 
TL model with a Peierls distortion also exhibits the prop- 
erties of a BI with equal charge and spin gaps. 



B. Renormalization effects from the short-range 
interaction 

In order to account for the renormalization of the BI 
gap by the short-range part of the interaction, we con- 
sider again the limiting case ?7 <C A <C t. Wc treat 
both the ionic term and the Hubbard interaction on the 
same footing in a pcrturbative way. For this purpose it 
is convenient to consider the 2D Euclidean action: 



(22) 



S — Sq + Sint , 

So=vfJ d^x {lliVcPcf + {Vcp^f]} , 
Sint= vf J d^x sin (^^/2^^ Kc<j)c^ cos (y2Tr<j)s^ 

+ ~l (a/SttK'c'^c) + cos (ySncpsj I • 
Here the dimensionless coupling constants are given by 
= M, 



U 



rur. 



Expanding the partition function 

Z = yD$«D$ee-^°[*-*=l[l-5i„t + i52„, + ...] (24) 

we integrate out all configurations in the second order 
term with |x — x'| ^ oq. Using the operator product 
expansion formulas 



/ ^sin(<^(x)).cos(2<^(x')) = 
^ / ^cos((/.(x)).cos(2(/.(x')) = 

TT f d^X / , / \x 1 

^2 J -^[cos(<;^c(x))-F...] , 

/ ^ / ^ COS (^(x)) • COS (^(x')) = 
J % J "0 



+ 



d^x. 



:(V(/)c)^ + cos(2^(x))-F. 



(25) 



(26) 



(27) 



A = A ( 1 - 



Mc = -U\^l- A, 
M, = u(l-X 



A2 
Ut 

Ut 



(29) 
(30) 
(31) 



A, Ac, and Ag are positive numbers of order unity. We 
observe that 

• the parameter of the ionic distortion is renormalized 
linearly by the on-site repulsion; 

• the amplitudes of the Umklapp and backward scat- 
tering processes (Mc and Mg, respectively) are renormal- 
ized quadratically by the ionic distortion; 

• the U —U asymmetry of the model is clearly man- 
ifested: at [/ > the amplitude of the strongly relevant 
Umklapp scattering decreases, while at [/ < the am- 
plitude of the strongly relevant backscattering processes 
increases. 

Therefore we conclude that in the IHM the BI phase 
is more stable for a repulsive than an attractive Hubbard 
interaction. 

For a Peierls distortion the effect of the electron- 
electron interaction at short distances is rather different. 
Assuming U <^ 6 <^t the same procedure as above gives 
in this case the following expressions for the renormalized 
parameters 



(32) 

(33) 
(34) 



Therefore, contrary to the ionic case, 

• the amplitude of the Peierls distortion is not renor- 
malized to linear order in U ; 

• the amplitudes of the Umklapp and backward scat- 
tering processes are renormalized quadratically by the 
Peierls distortion term; 

• at ?7 > the amplitude of the strongly relevant Umk- 
lapp scattering increases, while at /7 < the amplitude of 
the strongly relevant backscattering processes increases. 

Based on these arguments, one is led to predict that 
for the PHM the deviation from the BI behavior due to 
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the on-site coupling of arbitrary sign could be similar to 
that of the attractive IHM. 

Therefore we conclude that the BI phase in the Peierls- 
Hubbard model is less stable against a Peierls distortion 
which is not renormalized to linear order in U, and the 
amplitudes of the relevant scattering processes at 
U > and at J7 < always increase. The devia- 
tion from the BI behavior in the case of the attractive 
IHM and the PHM should therefore happen almost at 
the same values of the on-site interaction. The funda- 
mentally different behavior of the repulsive IHM and the 
PHM will indeed be verified below. 



C. The FGN phase diagram 

In Ref. 1^ FGN proposed a new scenario for the BI 
to MI crossover in the IHM. The key ingredient of their 
theory is the presence of two separate transitions: an 
Ising-type transition at U^/^ where the charge gap van- 
ishes and a continuous transition at U^p > C/f^ where the 
spin gap vanishes. The charge excitations are gapped for 
arbitrary U ^ J7^^, while the spin sector remains gapless 
at t/ > U^p. The ground-state phase diagram of the IHM 
and the properties of the charge and spin gapped phases 
were argued to be qualitatively captured by the effective 
potential: 

V "A sin (^^/2t:K,^c) cos (\/2^<?!>,) (35) 



Mc cos VSttK, 



AL cos 



where Mc,Ms,A > are phenomenological parame- 
ters obtained by integrating out high-energy excitations. 
Minimizing V (0c, 0s) with respect to </>c and 4>s results in 
the following sets of vacua: For A > 4Mc (i.e. U < C/^^) 



I. V27r(7 



— (mod 27r) 
2 ^ ' 



II. \/27r0s = TT + 2?™, ^2TxKc(i)c = - \ (™od 1-k) (36) 



while for A < 4Mc (i.e. U^^^ <V <U\ 



spl 



I. \/2tt4>s — 'i.Tin , 2ttKc4'c — VOi — V'o (mod 2tt) 
II. V2TT(f)s = TT + 27rn , 

\/2ttKc 4>c = — fo, — + fo (mod 27r) (37) 

where (po = arcsin(A/4Mc). 

The symmetry properties of the various ordered phases 
are described by the order parameters for the short wave- 
length fluctuations of the 

• siie-located charge density wave: 



AcDw =^ i-irY.' 



.t 



A 



BOW P 



CDW 



BOW 
CDW 



FIG. 1. Illustration of the CDW -h BOW phase. Top figure: 
two degenerate dimerization patterns; the vertical dashed line 
is a mirror axis indicating the link-reflection symmetry of each 
dimer state. Middle figure: schematic representation for a 
staggered CDW which has site-reflection symmetry. Lower 
flgure: the superposition of the CDW state with both dimer- 
ized patterns again allows for a site- reflection symmetric state. 



• siie-located spin density: 



- cos(\/27r/Cc0c) sm{V2^(l)s) , (39) 
• ftonrf-located charge density wave: 

Abow = (-l)"^(4^^c„+i,<T + h.c.) 



cos 



{^^271 Kc(j)c) cos{V2tt(I)s) 



(40) 



sin{\/2TrK^(/)c) cos(\/27r(/)s) 



(38) 



For U < U^i^ the charge and spin excitation spectra are 
gapped. The vacuum values of the ordered fields are 
{(j)s) = and (0c) = \/t^/&Kc and the system has long- 
range CDW order. Thus the set of vacua (|3|) corre- 
sponds to the BI phase. 

For U^i^ < U < Ugp the vacuum expectation values of 
the ordered charge field is now different: {(f) s) — and 
{\/2-kKc4>c) = Vo- In this phase, long-range CDW and 
BOW correlations coexist, 

{AcDw{x)AcDw{x')) ^ sin^ (po , 
{Abow{x)Abow{x')) --cos^<^o ■ 

The standard BOW phase violates site- but preserves 
link-inversion symmetry. In the CDW + BOW phase, 
due to the CDW pinning by the ionic distortion, the 
link-parity is also broken. However, in the absence of 
an externally imposed Peierls distortion, the BOW pat- 
tern is doubly degenerate. Thus the charge distribution 
in the CDW -I- BOW phase can be represented as a lin- 
ear combination of two dimerized patterns shifted with 
respect to each other by one lattice spacing (see Fig. |l|) 
for which the site-inversion symmetry is not broken. 
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FIG. 2. Lowest energy eigenvalues of the IHM at half-filling 
for L = 8 sites, periodic boundary conditions and A = O.bt. 



Finally, at J7 > U^p the spin gap disappears. According 
to FGN at this point (fo = and the standard MI phase, 
with an identical power-law decay of the SDW and Peierls 
correlations is realized. 

However, the results of our numerical studies presented 
below as well as QMC calculations and finite clus- 
ter studies using valence bond techniques jl^ indicate a 
different scenario. The amplitude of the CDW correla- 
tions smoothly decays with increasing U with ipQ for 
U ->oo. 



V. LANCZOS EXACT DIAGONALIZATION 
RESULTS 

In order to explore the nature of the spectrum and the 
phase transition, we have diagonalized numerically small 
systems by the Lanczos method |p3| extending earlier 
exact diagonalization calculations ||l0[,p| . The energies of 
the few lowest eigenstates were obtained for finite chains 
with L = 4n and PBC or L = An + 2 with APBC, for 
reasons discussed above. 

We first analyze a short chain without finite-size scal- 
ing. For chain lengths L < 16 finite-size effects do not 
change the qualitative behavior discussed below. In Fig. 
^, the lowest eigenenergies of the IHM for A = 0.5i, L — 8 
and PBC are shown as a function of U. At U = 1.3t, a 
level crossing of the two lowest eigenstates occurs. A non- 
degenerate eigenstate of the IHM has site-parity eigen- 
values ±1, so a ground state level-crossing transition cor- 
responds to a change of the site-parity eigenvalue. 

For U = 0, the IHM is easily diagonalized in momen- 
tum space by introducing fermionic creation (annihila- 
tion) operators 7^^^ ilkab) with an index b — 1,2 de- 
noting the lower and upper bands, respectively, giving 
two energy bands £^1/2 (fc) = ±^4cos^(fc) -I- (A/4)2 with 



FIG. 3. Overlap matrix elements of the exact ground states 
lipo) = \U/t,A/t,S} of the Peierls chain ([/ = A = 0) and 
the Peierls-Hubbard model (stars), the ionic band insulator 
(A = 0.5t, S = 0) and the IHM (full circles) and of the IHM at 
a given U and at U = 2.9t (circles) as a function of U. Calcu- 
lations were performed for 10 sites with antiperiodic boundary 
conditions. 



momenta — 7r/2 < k < Tr/2. For U = the first two 
degenerate excited states at half-filling always have neg- 
ative site parity, because the ground state has P = +1, 
and the operator 7^g.27g(Ti with q ~ 7r/2 obeys 



Pryt 

Inc 



q(72^qal 



7gcr27gcrl-P 



(41) 



The first two excited states shown in Fig. g are the spin 
singlet {S = 0, Sz — 0) and triplet excitations {S = 1, 
Sz — 0), created from the ground state by applying the 
operators 



^(7lT27,Ti+7^i27,xiy , 



(42) 



respectively. Thus both excited states have total momen- 
tum kfot = and negative site parity. For C/ > 0, these 
degenerate excited states split in energy. 

Fig. 1^ shows the overlap matrix elements between the 
exact ground states |V'o) = \U/t, A/t, 6) of the IHM, the 
PHM, and the Hubbard model. Specifically we show 
the overlap integrals between the ground states of the 
IHM for < f7 < 2.9t and U = 0, and between the 
ground states for < [/ < 2.9t and U = 2M. Clearly, 
a sharp transition occurs at Uc = 1.5t with the groimd 
states I [7,0.5,0) for U < Uc and U > Uc being orthog- 
onal to each other. This reflects the level-crossing phe- 
nomenon with a site-parity change. For the PHM (A = 0, 
d = 0.5) the situation is different. The overlap integral 
between the ground states for J7 = and U > decreases 
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continuously with increasing U, and no indication of a 
phase transition is observed. This model is a BI only at 
U ^ 0, and turns into a correlated Peierls insulator for 
any U > 0. 

Obviously, exact diagonalization of finite rings identi- 
fies one critical Uc > 0, separating a BI with P = +1 at 
U <Uc from a CI with F = -1 for [/ > Uc. 



VI. DMRG RESULTS 

In order to access the transition scenario in the long 
chain-length limit, we have studied chains up to L=300 
using the DMRG method ||,||. The fact that the 
transition at Uc is connected with inversion symmetry 
requires some caution when open boundary conditions 
(OBC) are used in our DMRG studies. For OBC and 
L = 2n the IHM is not reflection symmetric at any site. 
Thus, the ground state does not have a well defined site- 
parity, and the level-crossing transition is absent. To 
overcome this problem, one might try to use chains with 
OBC and an odd number of sites L = 2n + 1, since the 
Hamiltonian in this case is reflection-symmetric with re- 
spect to the site ic in the center of the chain, and a site- 
parity operator is well defined by 



Pc! „P^ 



-L+l- 



(43) 



To test whether this is an improved choice we have cal- 
culated the site-parity of the ground state for U — 
analytically for different chain lengths L = 2n + 1 and 
found 



P|Vo) = (-l)"|^o) 



(44) 



On the other hand, in the large U limit again the mapping 
to the Heisenberg Hamiltonian could be used to obtain 
the site-parity eigenvalue. By extending the idea of Gi- 
dopoulos et al. |^ to chains with L — 2n + 1, for U '3> t 
we obtain 



P\^o) = (-1) 



1) [EL'i H 



IVo) = (-l)"|V'o) 



(45) 



Thus, the parity eigenvalue of the ground state is the 
same at J7 = and U ^ t for a given chain length, and 
no level crossing occurs. As mentioned above, the ar- 
gument for J7 3> i requires some caution for the IHM, 
so we also have checked this result numerically by exact 
diagonalization of small systems with L = 5, 7, 9, . . . to 
confirm the absence of a level crossing. Obviously, results 
for any choice of boundary conditions should recover the 
level-crossing scenario when extrapolated to the thermo- 
dynamic limit. Due to the fact that the sharp transition 
at a well defined Uc does not exist in the finite-chain re- 
sults for OBC, the extrapolation is a rather subtle prob- 
lem, since a sharp transition feature has to be identified 
from the extrapolation of smooth curves. This requires 
the use of quite long chains in the critical region. 




FIG. 4. Results for the spin (As) and charge (Ac) 
gaps of the IHM at half-filling with A = 0.5t as a func- 
tion of U. Energies were obtained by DMRG calculations 
on open chains with L = {30, 40, 50, 60} (main plot) and 
L = {30, 40, 50, 60, 200, 300} (inset), and extrapolated to the 
limit of infinite chain length. 



We have evaluated Ac and A5 for chains up to 300 
sites. In Fig. || extrapolated results are shown as a func- 
tion of U. Calculations were performed with OBC for 
chains of lengths L = {30,40,50,60}, and additionally 
for L — 200 and L — 300 in the transition region around 
the estimated Uc- Unlike the definition (|^), Ac was ob- 
tained here using 



Ar 



1 



EoiN ^L + 2,S,^0) 



Eo{N = L-2,S, = 0)- 2Eo{N = L,S,=Q) . 



in order to calculate the relevant energies Eq in the sub- 
space Sz = 0. This becomes equivalent to (||) in the 
thermodynamic limit. We assume a scaling behavior of 
Ac and A5 of the form p& 



A,(L) 



L 



El 
L2 



(46) 



where i € {S,C'}. The extrapolation for L — > 00 is 
then performed by fitting this polynomial in 1/L to the 
calculated finite-chain results. We note that different 
finite-size scaling formulas were proposed in the litera- 
ture mainly when PBC or APBC were used 

As can be seen from the main plot in Fig. |, extrap- 
olating the results for L = {30,40,50,60} does indeed 
not give a sharp transition behavior of the gaps. As il- 
lustrated in the inset, adding results for L = 200 and 
300 in the critical region changes the picture consider- 
ably. Within numerical accuracy the charge and spin 
gaps remain equal up to a critical Uc~2.lt. A sharp 
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FIG. 5. Results for the spin (As) and charge (Ac) gaps of 
the Peierls-Hubbard model at half-filling with a modulation of 
the hopping amplitude 5 = 0.5 as a function of U. Energies 
were obtained by DMRG calculations on open chains with 
L = {30,40,50,60}. 

kink for Ac is observed reflecting the merging of the 
ground states of the different site-parity sectors. Impor- 
tantly, Ac does not close at the critical point. This is 
in fact not in conflict with an underlying ground-state 
level crossing. If the ground states of the different site- 
parity sectors become degenerate, the only rigorous con- 
sequence is the closing of the optical excitation gap. The 
selection rules for optical excitations (see section VII) 
allow only for transitions between states of different site- 
parity. Furthermore, optical transitions occur within the 
same particle number sector. The optical gap is there- 
fore by definition distinct from the charge gap Eq. (||) 
which involves the removal or the addition of a particle. 
The critical point Uc of the IHM corresponds to the re- 
markable situation where the optical gap closes while Ac 
remains finite. This implies that each site-parity sector 
separately has a finite Ac. Above Uc the charge and 
spin gaps split indicating that the corresponding insulat- 
ing phase is no longer a BI. As continuously decreases 
with increasing U and becomes unresolvably small within 
the achievable numerical accuracy. 

Our results are in general agreement with the data 
obtained by Qin et al. These authors performed 

DMRG calculations for the IHM with A = 0.6t, for 
chains up to i = 600 sites. In contrast to our calcula- 
tions, Qin et al. used the definition (|^) to calculate Ac. 
Surprisingly, they observed a non-monotonic scaling be- 
havior of As with L for values of U close to the critical 
t/c, i.e. for chain lengths L > 300 As started to increase 
again. It remains unclear whether this is due to loss of 
DMRG accuracy with increasing chain lengths. We note 
that in their substantially revised paper Qin et al. have 
presented additional DMRG data which agree with our 




FIG. 6. The electron density distribution in the ground 
state of the IHM for A = 0.5t and U = 0.8f (diamonds) and 
U = Ai (circles). Results were obtained by DMRG calcula- 
tions on open chains with L = 32. Stars show results for the 
Hubbard chain (A = 0). 

present finding. DMRG calculations for the IHM with 
A = t have also been performed by Takada and Kido for 
chains up to L = 400 sites The authors interpret 

their results in the region close to Uc in favour of a two- 
transition scenario similar to that of FGN |^ . Below we 
will show that such an interpretation is not valid. 

For comparison we show in Fig. ^ the spin and charge 
gaps versus U in the PHM. As already anticipated in pre- 
vious sections this model does not show any signature of 
a phase transition; i.e. Ac > As > for all J7 > 0. So 
although the Peierls and the ionic insulator for U — Q 
similarly possess an excitation gap at the BZ boundary, 
applying a Coulomb U leads to distinctly different be- 
havior in both cases. The origin of the different behavior 
must be traced to the fact that the Hubbard interac- 
tion and the ionic potential compete locally on each site, 
while the Peierls modulation of the hopping amplitude 
tends to move electronic charge to the bonds between 
sites, thereby avoiding conflict with the Hubbard term. 
This difference is also reflected in the characteristic struc- 
tures of the ionic and the Peierls terms in the bosonized 
versions of the Hamiltonian. 

The important question remains about the nature of 
the insulating phase of the IHM for U > Uc- As we 
argued above the numerical search for the vanishing of 
the spin gap will remain ambiguous due to the necessity 
to rule out the possibility of a finite but exponentially 
small As. To further analyze the BI and CI phases below 
and above Uc, we have calculated correlation functions 
in the ground state of the IHM using DMRG results for 
finite chains. 

Fig. H shows the charge distribution (0|(nr,| + n^^i — 
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U=0.8, A=0.0 E 




FIG. 7. Spin-spin correlation function in the ground state 
of the IHM for A = 0.5t (open symbols) and the Hubbard 
model (A = 0) (full symbols) at U — 0.8t (diamonds) and 
U ^4t (circles). Chain lenght L = 32. 




FIG. 8. Bond density correlation function in the ground 
state of the IHM for A = 0.5t (open symbols) and the Hub- 
bard model (full symbols) at [/ = 0.8t (diamonds) and U — it 
(circles). 



1)|0) in the ground state of the Hubbard model and the 
IHM at t/ = 0.8t < Uc and U = U > U, for a L = 32 
chain. The alternating pattern in the density distribution 
is well pronounced not only in the BI phase but also in the 
CI phase far beyond the critical point at J7 ^ Uc- For the 
L = 32 chain the CDW is well established at distances 
I ^ L/2 even at J7 = At. We note, that the attraction 
(or repulsion) of the charge from the chain edges, a well 
pronounced boundary effect in the BI phase, is absent at 
A = and in the CI phase of the IHM. The amplitude of 
the CDW pattern smoothly decreases with increasing U . 
Our numerical data indicate that the alternating pattern 
in the electron density distribution in the IHM remains 
for arbitrary finite U . 

Fig. shows the DMRG results for the spin-spin cor- 
relation function (0|S"=(L/2)S'^(L/2 + r)|0) in the ground 
state of the IHM at = 0.8t < Uc and U = U > Uc 
in comparison with the spin correlators in the Hubbard 
model. In the BI phase at U = 0.8t the SDW correla- 
tions are almost completely suppressed on a scale of half 
the chain length. At U — At the amplitude of the SDW 
correlation in the CI phase of the IHM is slightly reduced 
in comparison to the Hubbard model at the same value 
of U . However, the large distance behavior of the spin 
correlations in the CI phase and the MI phase of the Hub- 
bard model is similar. Therefore at arbitrary U > Uc the 
spin correlations in the CI phase of the IHM are almost 
identical to the Hubbard model. This equivalence, how- 
ever, is valid up to the accuracy of an unresolvably small 
gap in the spin excitation spectrum. 

To address the dimerization tendencies in the CI phase 
we have calculated the ground-state distribution of the 
normalized bond density 



BDir) = TTT - 1 • (4^) 

Fig. H shows the results of the DMRG calculations for a 
i = 32 IHM chain and the Hubbard chain at [/ = 0.8i 
and U = At. The boundary effect of an open chain is 
strong and leads to a modulation of the bond density al- 
ready for the pure Hubbard model. The comparison with 
the Hubbard chain indicates, that the ionic distortion 
leads to a weak suppression of the bond-density ordering 
at [/< C/c, while at C/ > Uc, the amplitude of the bond- 
density modulations slightly increases in the CI phase. 

To study the tendency in the IHM towards BOW or- 
dering at the transition into the CI phase we compare 
the boundary induced alternating patterns of the bond 
density at A = and A = 0.5t as a function of U. For 
the open L = 32 IHM chain the effect is clearly seen if we 
compare the corresponding bond densities in the central 
bond of the chain as shown in Fig. ||. In the BI phase at 
U < Uc, the ionic term reduces the bond density. How- 
ever in the vicinity of the transition point and at U > Uc 
the bond density in the IHM is larger than in the pure 
Hubbard model. It is notable that this difference remains 
positive for arbitrary U > Uc and disappears only in the 
limit U oo. 

We conclude that the finite chain DMRG studies of the 
IHM indicate the presence of only one transition from 
the BI to the unconventional CI phase with long range 
CDW order for all [/ < oo. In the CI phase the spin 
excitation spectrum is either gapless or characterized by 
an unresolvably small spin gap. If the spin gap is finite, 
we expect true long range BOW order coexisting with a 
CDW. If however the spin gap vanishes, SDW and the 
dimer-dimer correlations coexist in this phase with an 
identical decay at large distances in the presence of a 
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FIG. 9. Difference between tlie bond-cliarge density of the 
IHM witli A = 0.5t and tfie Ifubbard model for A = as a 
function U. Tc denotes the bond {rc,rc + 1) in the center of 
a L = 32 chain. 

finite CDW amplitude. The fact, however, that the spin 
gap is finite and equal to at the transition supports 
the scenario of a CI phase for U > Uc with a finite As 
and coexisting long range CDW and BOW order. 



VII. OPTICAL CONDUCTIVITY 

To further analyze the BI and CI phases we have cal- 
culated the frequency-dependent conductivity within the 
dynamical DMRG (DDMRG) approach, making use of 
the correction- vector technique described in Rcf. 
Calculations were performed for L — 128 sites with OBC 
using a Parzen filter to suppress the influence of the chain 
boundaries. 

The real part of the dynamical conductivity cr(w) is 
determined from the linear response of the system to an 
external electromagnetic field. In the Kubo formalism, 
cr{uj) is related to the imaginary part of the retarded 
current-current correlation function 



La 



dre^'^^ (^o|[j-,(r),j,(0)]|^o) 



E 



\{^n\jg\i^0)\' 



La ^^\hLO + {En - Eo) + iO+ 

l(^n|j,|V'0)P 

;iw-(K--Bo) + iO- 



with T the time index, a the distance between neighbor- 
ing sites, and |V'n) and E^ denoting the eigenstates and 
their respective energies. The paramagnetic current op- 
erator is defined by 



3q 



h ^ 

La 



(49) 



where e is the electron charge. The structure of the cur- 
rent operator and the required matrix elements enforces 
the important selection rule that only transitions between 
states with different site-parity are allowed. The real part 
of the conductivity in the long-wavelength limit g = is 
given by 



reg 



> 0) = ImXjj (g = 0, u) 



(50) 
(51) 



For the insulating phases of the IHM at half-filling the 
Drude weight either vanishes, D = 0, or is not defined 
for degenerate ground states, so that cri(w) — a\'^^{uj). 

In the absence of interactions J7 = 0, the conductiv- 
ity diverges as <Ji(uj) ^ ^l\J^ — Ap for w Ac and 
a; > Ac We expect this behavior to persist upon 

increasing U over the entire BI phase. In order to test this 
expectation we have evaluated the optical conductivity in 
the BI phase by DDMRG . For the correction- vector 
method a separate DMRG run has to be performed for 
each selected frequency uj. In Fig. |l^, the DDMRG re- 
sults for a chain with L — 128 sites and OBC are plotted 
for \J = Zt. For the proper convergence of the correction- 
vector algorithm a finite broadening r/ = Q.lt was chosen 



in Eq 
between t 



\ line If 

- , (48) ^ 



For A = the solid black line interpolating 
die discrete data points at different frequencies 
was calculated by additionally applying the Lanczos vec- 
tor method after the DMRG sweeps, obtaining the corre- 
lation function also in the vicinity of the selected frequen- 
cies . One observes a dominant low-energy excitation 
peak above the optical gap at around E « l.lt. The 
Lorentzian tail for E < 1.7t results from the necessarily 
large broadening. Also, in the correction-vector results 
the square-root divergence at the absorption edge is not 
visible due to the finite size of the system. However, it 
is to be anticipated from the dominant excitation peak 
above the gap edge. 

For U > Uc, the situation is different, as demonstrated 
in Fig. |l^, too, where we also show (Ti(w) for U = 3t 
and A = Q.5t. The dashed line connecting the data 
points serves only as a guide to the eye. In the CI phase, 
the behavior of the optical excitations is expected to be 
similar to those of the repulsive Hubbard model, where 
ai{uj) ~ — Ac above the charge gap Ac which in 
this case is identical to the optical gap |Q. In the in- 
set of Fig. |o|, (Ti{uj) is plotted on a finer energy scale 
in the onset region above the gap edge. The solid black 
line is a fit to the DDMRG data obtained by convoluting 
with a Lorentzian of width 77, i.e. by using 



^opt 



the fit formula 



A 



^opt 



(52) 



where A is an adjustable prefactor and Vmax the upper 
bound for the square-root dependence. The fit in Fig. |l^ 
has been obtained using the parameter values A — 1.9t, 
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FIG. 10. Real part of the optical conductivity cri(tj) in 
units of Roa = e^a/h vs. energy E = hco for the half-filled 
IHM for U = 3i. The results were obtained by DDMRG 
calculations using the correction-vector method for a chain 
with L = 128 sites and OBC. A finite broadening rj = O.lt 
has been used. The selected A-values correspond to the BI 
(A — 4t) and the CI (A = 0.5t) phase. The inset shows a fit 
(solid black line) to the onset region of ai{Lj) in the CI phase. 

= l.bt, and Aopt = 0.48i. The value for Aoct 
is close to the charge gap result presented in Fig. 
A difference between Aopt and Ac can not be resolved 
within the accuracy of the DDMRG data. The agree- 
ment between the fit and the DDMRG results is quite 
striking. Thus, we conclude that the numerical results 
are in agreement with the behavior expected for a cor- 
related (Mott) insulator in ID. The qualitative behavior 
of cri(w) compares well with the results obtained for the 
Hubbard model at U = 3t 

Further consequences for the optical conductivity fol- 
low from the level crossing scenario discussed in the previ- 
ous sections. On approaching Uc the ground and the first 
excited state become degenerate, Ei —^ Eq; the matrix 
element {i^i\jq=o\4'o) between these states remains finite 
due to their different site parity eigenvalues. Therefore, 
the optical gap vanishes precisely at Uc- The critical 
point has been commonly referred to in the literature as 
"metallic" . However, by definition a metal is character- 
ized not only by gapless excitations, but also by a finite 
Drude weight D > 0. Due to the ground-state degener- 
acy the Drude weight for U — Ucis not properly defined. 
Thus it follows that although the optical excitation gap 
vanishes in the IHM for U = Uc, it cannot be considered 
"metallic" in the usual sense. (Following the definition 
of W. Kohn a metal has a finite Drude weight, D > 

E3-) . 

We emphasize that the optical conductivity was so far 
evaluated in the BI and the CI phase sufficiently far from 
the critical point. Our currently achievable accuracy does 



not allow safe conclusions for the optical conductivity in 
the critical region around Uc- At Uc the charge gap Ac 
remains finite while the optical gap vanishes. This spe- 
cial situation demands that near below Uc the optical 
conductivity has an isolated peak below the onset of the 
continuum for optical absorption. The energy of the iso- 
lated peak corresponds to the energy difference of the 
individual ground states in each site-parity sector with 
P = ±1. Resolving this remarkable feature will remain 
a challenge for future numerical simulations. 

VIII. CONCLUSIONS 

Our combined analytical and numerical analysis clari- 
fies the ground-state structure of the IHM. A single tran- 
sition at a critical /7c (A) separates a CDW-BI from a 
CI phase with finite charge and optical excitation gaps 
and most likely also a finite spin excitation gap. Above 
Uc, CDW and BOW order coexist and persist for all 
U > Uc- No secondary transition into a true MI phase 
occurs which - in fact - would have a higher symmetry 
than the IHM itself. The insulator-insulator transition 
on finite chains of the IHM results from a ground-state 
level crossing of the site-parity sectors with P = ±1. The 
BI phase has a unique ground state with P — +1 while 
the ground state in the CDW -I- BOW phase above Uc is 
necessarily two- fold degenerate. The critical point Uc{A) 
is characterized by the remarkable situation that the op- 
tical absorption gap vanishes while the spin and charge 
excitation gaps remain finite and equal. The distinction 
between the optical and the charge gap is the key feature 
for the structure of the insulating phases of the IHM. Re- 
fined numerical studies of the optical conductivity in the 
critical region define therefore the future task. 

The existence of an insulator-insulator transition in 
the IHM at half-filling naturally raises the question of 
whether a ground state phase transition finds its contin- 
uation also when the system is doped away from half- 
filling. The possible consequences for the dynamics of 
doped charge carriers and for the phenomenon of spin- 
charge separation will be the topic of future investiga- 
tions in the IHM. 
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